% get_duration_loadwritewrite
% 采集每天指定时间段内的负荷
load("yantian402\2024年配变时序数据.mat","result2");
all_load_data=result2;
clear result2

end_date=max(all_load_data.Date);
begin_date=min(all_load_data.Date);
begin_time_per_day=duration(0,0,0);
end_time_per_day=duration(23,59,59);
zero_clock=dateshift(begin_date,'start','day');
begin_ticks=zero_clock+begin_time_per_day:duration(24,0,0):end_date;
end_ticks=zero_clock+end_time_per_day:duration(24,0,0):end_date+duration(24,0,0);
Nt=length(begin_ticks);

% 分箱
Y1 = discretize(all_load_data.Date,[begin_ticks end_ticks(end)]);
Y2 = discretize(all_load_data.Date,[begin_ticks(1) end_ticks]);
idx=(Y1==Y2);
all_load_data=all_load_data(idx,:);
all_load_data.timebins=Y1(idx);
[duration_P,avg_P,max_P,min_P,std_div_P,cnt]=splitapply(@extract_info,all_load_data.P,all_load_data.timebins);

NsampleDays=15;

stat_table=table(duration_P,avg_P,max_P,min_P,std_div_P,cnt);
stat_table.timebins=reshape(1:Nt,[],1);
%过滤stat_table，仅留下cnt在[min(cnt)+0.1*(max(cnt)-min(cnt)),max(cnt)]区间内的数据,算临时修改
% valid_range = [min(stat_table.cnt) + 0.1 * (max(stat_table.cnt) - min(stat_table.cnt)), max(stat_table.cnt)];
% stat_table = stat_table(stat_table.cnt >= valid_range(1) & stat_table.cnt <= valid_range(2), :);

%将每天平均时点负荷分箱
% histogram(avg_P,NsampleDays);
[~, avg_P_edges,avg_P_bins] = histcounts(stat_table.avg_P, NsampleDays);
stat_table.avg_P_bins=avg_P_bins;

%将每天的极差分箱
% histogram(max_P-min_P,NsampleDays);
[~, max_min_P_edges,max_min_P_bins] = histcounts(stat_table.max_P - stat_table.min_P, NsampleDays);
stat_table.diff_P_bins=max_min_P_bins;

%将每天的标准差分箱
% histogram(std_div_P,NsampleDays);
[~, std_div_P_edges,std_div_P_bins] = histcounts(stat_table.std_div_P, NsampleDays);
stat_table.std_div_P_bins=std_div_P_bins;


%手动挑十个，画出日负荷阶梯图
select_table=struct('desc',["0410","0203","0728","1029","夏至0621","除夕0209","国庆节1001"].',...
    'indices',[101,34,210,303,163,30,276].');
select_table=struct2table(select_table);


% 绘制年负荷阶梯图
figure('Units','centimeters','Position',[5,9,20,7]);
stairs(1:height(stat_table),stat_table.duration_P, 'LineWidth', 1.5);
xlabel('天数');
ylabel('总负荷 (kW)');
xlim([1,366]);
title('2024年负荷曲线');

% 定义2024年四季的分界点（以立春、立夏、立秋、立冬为准）
season_dates = {'2024-02-04', '2024-05-05', '2024-08-07', '2024-11-07'}; % 立春、立夏、立秋、立冬
season_names = {'立春', '立夏', '立秋', '立冬'};
season_colors = {'g', 'r', [0.9290 0.6940 0.1250], 'b'}; % 绿、红、橙、蓝代表四季
% 计算分界点对应的天数
base_date = datetime('2024-01-01');
season_days = zeros(1,4);
for i = 1:4
    current_date = datetime(season_dates{i});
    season_days(i) = days(current_date - base_date) + 1;
end
% 标记四季分界点
hold on;
for i = 1:4
    x = season_days(i);
    plot([x x], ylim, '--', 'Color', season_colors{i}, 'LineWidth', 1);
    text(x, max(ylim)*0.95, season_names{i}, ...
        'Color', season_colors{i}, ...
        'HorizontalAlignment', 'center', ...
        'BackgroundColor', 'w');
end
% 标记季节区域
text(mean([1 season_days(1)]), max(ylim)*0.85, '冬季', ...
    'Color', 'b', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
text(mean([season_days(1) season_days(2)]), max(ylim)*0.85, '春季', ...
    'Color', 'g', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
text(mean([season_days(2) season_days(3)]), max(ylim)*0.85, '夏季', ...
    'Color', 'r', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
text(mean([season_days(3) season_days(4)]), max(ylim)*0.85, '秋季', ...
    'Color', [0.9290 0.6940 0.1250], 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
text(mean([season_days(4) 366]), max(ylim)*0.85, '冬季', ...
    'Color', 'b', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');

% 标记重要节日 - 使用圆点标记
% 1. 标记除夕 (2024年2月9日)
chuxi_date = datetime('2024-02-09');
chuxi_day = days(chuxi_date - base_date) + 1;
chuxi_load = stat_table.duration_P(chuxi_day);
plot(chuxi_day, chuxi_load, 'ro', 'MarkerSize', 8, 'LineWidth', 1.5, 'MarkerFaceColor', 'r');
text(chuxi_day, chuxi_load*0.9, ' 除夕', 'Color', 'r', 'FontWeight', 'bold');
% 2. 标记国庆节 (2024年10月1日)
national_date = datetime('2024-10-01');
national_day = days(national_date - base_date) + 1;
national_load = stat_table.duration_P(national_day);
plot(national_day, national_load, 'bo', 'MarkerSize', 8, 'LineWidth', 1.5, 'MarkerFaceColor', 'b');
text(national_day, national_load*1.1, ' 国庆', 'Color', 'b', 'FontWeight', 'bold');

hold off;



% 绘制日负荷阶梯图
n_samples=height(select_table);

for i = 1:n_samples
    figure('Units','centimeters','Position',[6+0.2*i,10+0.2*i,16,7]);
    % 获取对应日期的数据
    day_index = select_table.indices(i);

    % 提取该天的数据,按照时点数相加
    day_data = all_load_data(all_load_data.timebins == day_index, :);
    [Grp,day_data_ticks]=findgroups(day_data.Date);
    day_data_P=splitapply(@sum,day_data.P,Grp);

    % 绘制阶梯图
    stairs(day_data_ticks,day_data_P, 'LineWidth', 1.5);
    xlabel('时间 (小时)');
    ylabel('总负荷 (kW)');
    title(select_table.desc(i));
    datetick('x', 'HH:MM', 'keepticks');
end

function [sum_d,avg_d,max_d,min_d,std_div,cnt]=extract_info(X)
    %应该按时间间隔算能量的

    cnt=numel(X);
    sum_d=sum(X);
    avg_d=mean(X);
    max_d=max(X);
    min_d=min(X);
    std_div=std(X);
end


